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The integral relations formalism introduced in [l|, 13 > an d designed to describe 1+ JV reactions, is 
extended here to collision energies above the threshold for the target breakup. These two relations 
are completely general, and in this work they are used together with the adiabatic expansion method 
for the description of 1+2 reactions. The neutron-deuteron breakup, for which benchmark calcula- 
tions are available, is taken as a test of the method. The s-wave collision between the 4 He atom and 
4 He2 dimer above the breakup threshold and the possibility of using soft-core two-body potentials 
plus a short-range three-body force will be investigated. Comparisons to previous calculations for 
the three-body recombination and collision dissociation rates will be shown. 

PACS numbers: 03.65.Nk, 21.45.-v,31.15.xj,34.50.-s 
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I. INTRODUCTION 

Calculation of continuum states corresponding to pro- 
cesses where a particle hits a bound iV-body system 
requires in principle knowledge of the corresponding 
(1 + jV)-body wave function at large distances. Need- 
less to say, the technical difficulties one has to face in 
order to obtain the wave function increase dramatically 
with N. In fact, already for N — 2, calculation of the 
three-body wave function is far of being trivial. 

However, even if knowledge of the (1 + iV)-body wave 
function is unavoidable, the possibility of reducing the 
distance at which such wave function is needed is in it- 
self an important step forward in the description of the 
reaction. As shown in [HQ, this can indeed be done 
by means of two integral relations that are based on the 
Kohn Variational Principle. These two integral relations 
are a generalization to more than two particles of the 
integral relation given in 0, 0] , and they permit to ob- 
tain the K,- (or <S)-matrix of the reaction by using only 
the internal part of the wave function. Therefore, all the 
physical information concerning a given 1 + N reaction 
can be obtained without an accurate knowledge of the 
asymptotic part of the wave function. 

The fact that the asymptotic part is not needed any- 
more leads to a drastic reduction of the computer effort 
required to extract the /C-matrix. The only condition 
necessary to obtain accurate second-order estimates of 
the /C-matrix through the integral relations is that the 
trial wave function must fulfill the Schrodinger equation 
only in the interaction region. This means that, for in- 
stance, scattering states can be described with bound- 
state-like trial wave functions [j| • 

In [l|, Q the integral relations were implemented to 
describe 1 + 2 reactions below the breakup threshold. 



When used in combination with the hyperspherical adi- 
abatic expansion method the dimension of the /C-matrix 
describing the reaction is dictated by the number of adia- 
batic potentials related to the possible outgoing elastic or 
inelastic channels, which typically is very small. Further- 
more, thanks to integral relations, the convergence of the 
/C-matrix in terms of the adiabatic channels included in 
the expansion of the wave function is highly accelerated. 
Actually, the pattern of convergence is similar to the one 
found when the same method is used for the description 
of bound states This is in fact not a minor issue, 

since as proved in [fj, when used to describe low-energy 
scattering states, the convergence of the adiabatic expan- 
sion slows down significantly, even to the point that an 
accurate calculation of the /C-matrix would require in- 
finitely many adiabatic terms in the expansion, what in 
practice makes the procedure useless. 

The success of the method for energies below the 
breakup threshold immediately suggests its extension to 
describe low-energy breakup reactions. In this case the 
dimension of the /C-matrix is not finite, since contrary 
to the elastic and inelastic channels or transfer reactions, 
the full three-body continuum states are described by in- 
finitely many adiabatic potentials. 

The first goal of this work is to generalize the method 
described in 

0,1 

for 1 + N reactions to energies above 
the threshold for breakup of the bound target. This gen- 
eralization will be shown in Section IT where a short 
summary of the method described in [2j will be given. 
In Section III, the n — d reaction, for which a series of 
benchmark calculations are available 0, [|| , will be used 
as a test of the method. 

The second goal, which will be discussed in Section 
IV, concerns the use of the new method to investigate 
the 4 He+ 4 He2 atomic reaction above the dimer breakup 



threshold. The case of the + state will be considered. In 
particular, we will focus on the use of soft-core 4 He- 4 He 
potentials. As shown in Q, the use of an attractive gaus- 
sian potential reproducing the same two-body properties 
as a standard hard-core potential (like for instance the 
LM2M2 potential) leads to equivalent bound three-body 
systems and phase-shifts for the elastic 4 He+ 4 He2 reac- 
tion, but only once that the soft-core potential is used 
together with a three-body short-term force. The possi- 
bility of using the same kind of potentials also to describe 
the breakup channel is interesting, since it automatically 
eliminates all the technical difficulties arising from the 
presence of a hard-core repulsion in the potential (see for 
instance Ref . [Tol| ) . 

Finally, we close this work with a short summary and 
the conclusions. 



II. FORMALISM 

In Refs.fl], 1 + 2 reactions were described within 
the framework of the hyperspherical adiabatic expansion 
method for energies below the threshold for breakup of 
the dimer. Therefore, only elastic, inelastic, and trans- 
fer processes were possible. Since the formalism is de- 
scribed in great detail in Ref. |2| , here we just summarize 
its main aspects, which are given in the first part of this 
section. In particular, we will focus on those key points 
that will permit an easier understanding of the gener- 
alization of the method to energies above the two-body 
breakup threshold, which will be shown in the second 
part of the section. In the last part we will describe 
the integral relations, which are actually the tools that 
permit to extract the /C-matrix of the reaction from the 
internal part of the wave functions. 



A. Sketch of the method for energies below the 
breakup threshold 

Given a three-body system, the corresponding wave 
function within the frame of the adiabatic expansion 
method is written as: 



P 1 



(1) 



where x and y are the usual Jacobi coordinates, and 
{p, £2} are the hyperradius and the five hyperangles ob- 
tained from x and y [ll|. The wave function has a well 
defined total angular momentum J, but for simplicity in 
the notation we omit this index (and its projection) in 
the expression above. 

In hyperspherical coordinates the Hamiltonian opera- 
tor H takes the form: 



h 2 f 

2m p 



(2) 
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FIG. 1: Typical effective adiabatic potentials for a three- 
body system where two two-body bound states are present. 
The two lowest adiabatic potentials go asymptotically to the 
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and E 2b of the two-body bound states. 
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binding energies E 2b 

For a given three-body energy E, when E^ < E < E 2b ' only 

(21 

one channel is open, when E 2b < E < both channels are 
open, and for E > the breakup channel is also open. 



where T p is the hyperradial kinetic energy operator, 7-Lq 
contains all the dependence on the hyperangles, and m 
is an arbitrary normalization mass. 

The angular functions $„(p, il), which form the com- 
plete basis used for the wave function expansion (fT]), are 
actually the eigenfunctions of Hq, 



in such a way that the adiabatic effective potentials 
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(4) 



enter in the coupled set of radial equations that permit 
to obtain the radial functions f n (p), 
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dp 2 
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where 



.(p) = (* n ( P) n)|JL|* n /(p,fi)) 



Pnn'(p) = ($ n (p,fi)|g-|*n'0>,n)) n , 



(5) 

(6) 
(7) 



where the subscript indicates integration over the hy- 
perangles only (see [ll| for details). 

A typical behavior of the adiabatic potentials is shown 
in FigfTJ They correspond to a three-body system where 
two of the two-body subsystems have a bound state. This 
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is reflected in the fact that the two lowest effective adia- 
batic potentials go asymptotically to the binding energies 
.E^b and E^ of each bound two-body system. 

In the figure the different regions defined by the energy 
of the incident particles are depicted. All the three-body 
energies E such that E^ < E < E$ (like E^ in the 
figure) correspond to processes where only one channel 
is open. Only the elastic collision between the third par- 
ticle and the bound two-body state with energy E^jJ is 
possible. When the three-body energy increases up to 
the region E$ < E < {E^ i n the fi gure) a second 
channel is open. Two different collisions are now possi- 
ble, the one where a particle hits the bound state with 
binding energy E$ , and the one where a particle hits 

(2) 

the state with binding energy EL ■ In the same way, 
each of these reactions has two possible outgoing chan- 
nels, corresponding to the two allowed bound two-body 
states and the third particle in the continuum. In other 
words, in this energy range the inelastic (if the two bound 
two-body states correspond to the same subsystem) or 
transfer (if the two bound two-body states correspond to 
different subsystems) channel is open. When E > (like 
E^ in the figure), the breakup channel is also open, and 
it is described by the remaining infinitely many adiabatic 
potentials. 

The coupled system of radial equations J5|) decouple 
asymptotically, and each of the radial wave functions be- 
have at large distances as dictated by: 



function is given by: 



2m dp 2 



e Up) = o. 



(8) 



When n corresponds to a closed channel the radial 
wave function /„ vanishes asymptotically. For values of 
E below the two-body breakup threshold, as considered 
in , and for a given incoming 1 + 2 channel, only out- 
going 1 + 2 channels are allowed (either elastic, inelastic, 
or transfer). In this case, as shown in (Tlj . the equation 
above describing the asymptotic behavior of the open- 
channel wave functions becomes 



Vn 



Up) = o, (9) 



where £ y is the relative angular momentum between the 
dimer and the third particle, y n refers to the modulus of 
the Jacobi coordinate between the center of mass of the 
outgoing bound two-body system and the third particle, 
and 



W = \lw< E - R - 



0)\ 
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(10) 



with E^' being the binding energy of the bound two- 
body system associated to the open channel n. 

With this in hand, it is not difficult to see that the 
asymptotic form of the corresponding three-body wave 



"0 



n=l 



(11) 



where i refers to the incoming channel, no is the number 
of open channels (all of them 1 + 2 channels), and 



T * n (p,n) (12) 



Gi K} = ^/kf } m y (k^y n )^n(p,n), (13) 



( ,3/2 

where jg, and rje are the usual regular and irregular Bessel 
functions (provided that we are dealing with short range 
potentials). 

As shown in [2], Eq. (|llj) can be written in a compact 
matrix form as: 



A( k )f (k) + B {K) G {K) = A {K) (f [k) - KG {K] 



(14) 

where \& is a column vector with no terms corresponding 
to the no open channels, A^ K > and B^ K > are no x no 
matrices made by the A^' and B\*p elements in Eq. (ITT1) , 
and F^ K > and G^ K > are again column vectors whose no 
terms are given by Eqs. (fT2"l) and (fl~3|). 

From the equation above it is clear that the /C-matrix 
of the reaction is given by: 



k = -aw~ 1 bW 



(15) 



This matrix is real, and from it one can easily obtain the 
S- matrix as (1 + i/C)(l — ilC)^ 1 . 



B. Generalization to energies above the breakup 
threshold 



When the total three-body energy E in a 1 + 2 reaction 
is above the threshold for breakup of the dimer, the first 
consequence is that infinitely many adiabatic channels 
are then open (see Fig. Q] for E — E^). Of course, still 
a finite number of them correspond to elastic, inelastic, 
or transfer processes, and the infinitely many remaining 
ones describe the breakup channel. Therefore, in this 
case the corresponding KL- (or S-) matrix has infinite di- 
mension. 

In any case, for the finite open channels corresponding 
to outgoing 1 + 2 structures, the expressions given in the 
previous subsection are still valid. This means that for 
these particular outgoing channels the three-body wave 
function behaves asymptotically as given by Eq. (fTTj) , and 
Eas.(fT2"]l and CE3]) are still valid. 

On the other hand, the breakup channels are char- 
acterized by the fact that the effective potentials V^J 
associated to them go asymptotically to zero as: 



ft 2 (*+§) (g 

2m p 2 



(16) 
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where K is the grand- angular quantum number defined 
as 2v + £ x + £ y , where £ x and £ y are the orbital angular 
momenta associated to the Jacobi coordinates x and y, 
respectively, and v = 0, 1,2, Therefore, asymptot- 
ically, each breakup adiabatic potential is associated to 
a fixed value of K. In fact, the corresponding angular 
eigenfunction $„(p, f2) is, also asymptotically, a linear 
combination of hyperspherical harmonics with that par- 
ticular value of K. 

When inserting p^|) into ([5]), we easily obtain that, 
asymptotically, the radial wave function for an outgoing 
breakup channel n satisfies the equation: 



dp 2 



(#+f )(#+§) 



Up) = o, 



(17) 



which is formally identical to the Eq.©, which is satisfied 
by the radial wave functions associated to outgoing 1 + 
2 channels, but replacing £ y by K + 3/2, k y n ^ by k = 
y / 2mE/h 2 , and y n by p. 

Therefore, in the more general case where the breakup 
channel is open, and assuming an incoming 1 + 2 channel 
i, the asymptotic form of the corresponding three-body 
wave function is given by: 



where F^' 



oo 

E 

n=l 



(18) 



and Gi^ 1 are still given by Eqs. (IT2l and 
(U2J), but of course, with the understanding that when n 
corresponds to an outgoing breakup channel the replace- 
ments t y — > K + 3/2, ky — >• K, and y n — >■ p have to be 
made (note that the relations j K+ a(z) = y^-j Jk-+2 (2) 
and i] K+ 3(z) — ^/j^Yk +2 (z) permit to write Eqs. (TT2")) 
and (|13[) in terms of the Bessel functions Jk+2 and Yk+2, 
which is how the asymptotic form of the breakup chan- 
nels is usually presented in the literature). 

Of course, the matrix form in Eq.(|14j) of the asymptotic 
wave function can still be used. As before, the /C-matrix 

of the reaction is given by K. = — B^ K \ but now 

the matrices A^ K > and B( K > have in principle infinite 
dimension and some truncation is then required. 

Due to the fact that the hyperradius p appears as the 
natural radial coordinate describing the asymptotic be- 
havior of the breakup channels, we then have that, for 
these channels, the adiabatic expansion is able to reach 
the correct asymptotic behavior for a sufficiently large, 
but finite, value of p. For this reason, for processes with- 
out 1 + 2 open channels (3—^-3 processes), the K.- matrix 
could in principle be extracted with sufficient accuracy 
from the asymptotic behavior. However, this is not true 
when 1 + 2 channels are open. For these channels the cor- 
rect asymptotic form (described by (12) and (13)) can not 
be reached until the Jacobi coordinate y and hyperradius 
p are equal, which only happens at infinity. Therefore it 
is essential in this case to use a formalism in which the 
/C-matrix is not extracted from the asymptotic part of 
the wave function but from its internal part. 



C. Integral relations 

The derivation of the integral relations has been shown 
in Ref.|2]. This derivation is completely general, and it 
is not particularized to the case of incident energies be- 
low the breakup threshold. Therefore, the same expres- 
sions derived in Q apply when the breakup channel is 
open. According to it, by making use of the Kohn Vari- 
ational Principle, it has been proved that when using a 
trial three-body wave function vp*, one can obtain the 
matrices A^ K > and B^ K ' accurate up to second order in 
(5(^ — \E''), and their matrix elements are given by: 



Bi 



2m 

2m 

'J 2 



^\n-E\F^ K) ) 



m\u-E\G 



(19) 
(20) 



where ^\ describes each possible incoming channel and 
the index j refers to each possible outgoing channel (ei- 
ther 1 + 2 or breakup) . 

In this work, the trial three-body wave function will 
be the one obtained as sketched in subsection III Al To be 
precise, it will be obtained by solving the Faddeev equa- 
tions by means of the hyperspherical adiabatic expansion 
method (see Ref.fi~l| for details). 

Since the regular and irregular functions F^ K ^ 
and G^ K } are asymptotically solutions of (H — 
E)F^ K \G {K) = 0, it is then clear that the integral rela- 
tions (HHJ) and ([2H]) depend only on the short-range struc- 
ture of the scattering wave function 

It is important to note that the function G^ K > , defined 
in Eq. (|13|) . is irregular at the origin. In order to avoid 
the problems arising from this fact, the G^ function is 
regularized. This means that in Eq. (|2"0| the Bessel func- 
tion rj contained in G^ is actually replaced by another 
function that goes to zero at the origin and behaves ex- 
actly as 77 at large distances. In Q this was done by using 
fjt(z) — (l — e~ lz y +1 rn(z), where 7 (> 0) is a parameter. 
However, in our case, where the index £ can reach pretty 
high values [l = K + 3/2 in the breakup channels) this 
procedure is not appropriate. 

In this work we have regularized the irregular Bessel 
function by solving Eq.® for each individual adiabatic 
potential. The solutions of this equation behave asymp- 
totically as fn(z) — > zj v (z) — tan S„ zr} v {z), where S n is 
the phase-shift, and the index v is equal to £ y for the 
1 + 2 channels, and K + 3/2 for the breakup channels. 
Then, the irregular Bessel function implicitly contained 
in Eq. (|2"0)) is replaced by: 



Zjv(z) - f n (z) 

z tan <L 



(21) 



which by construction is regular at the origin and goes 
asymptotically to r\ v (z). 
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FIG. 2: (Color online) Effective adiabatic potentials for the 
neutron-neutron-proton system in the doublet case. The thick 
curve corresponds to the neutron+deuteron channel, and it 
goes asymptotically to the deuteron binding energy. All the 
remaining adiabatic potentials correspond to breakup chan- 
nels. 



III. A TEST CASE: NEUTRON-DEUTERON 
SCATTERING 

The benchmark solutions for neutron-deuteron 
breakup amplitudes are shown in 0, Q • For this reason 
we take this case as a test for the method shown in this 
work. 

In particular, the nucleon-nucleon interaction is chosen 
to be the revised Malfliet-Tjon I-III model s-wave poten- 
tial Q, which for the spin triplet and singlet cases takes 
the form: 

V t (r) = - (-626.885e~ 1 - 55r + 1438.72 e - 311r ) (22) 
V s (r) = - (-513.968e- 1 ' 55r + 1438.72 e - 311r ) ,(23) 

where r is given in fm, and the potential in McV. Also, 
H 2 /m=41A7 MeV fm 2 . The potential above leads to a 
binding energy for the deuteron of 2.2307 MeV. 

In the calculation only s-waves are considered. There- 
fore, two different total angular momenta are possible, 
the quartet case (J = 3/2), for which only the triplet 
s-wave potential (1221) enters, and the doublet case (J = 
1/2), for which both, the singlet and the triplet potentials 
contribute. In Fig|2]we show the computed adiabatic ef- 
fective potentials (Eq.@) for the doublet case. As you 
can see, they follow the general trend of the potentials 
shown in FigfTJ although in this case there is only one 
1 + 2 channel, which corresponds to the neutron-deuteron 
reaction that we want to investigate. The corresponding 
adiabatic potential is given by the thick-solid line in the 
figure, and it will be labeled as channel 1. As expected, 
the asymptotic value of this potential corresponds to the 
binding energy of the deuteron. 



TABLE I: Inelasticity parameter |<Sn| for the neutron- 
deuteron scattering for two different laboratory neutron beam 
energies (14.1 MeV and 42.0 MeV) for the doublet and quar- 
tet cases. The value of K maK is the K-value associated to 
the last adiabatic potential included in the calculation. The 
numbers within parenthesis have been obtained without use 
of the integral relations. The last row gives the value quoted 
in Ref.[|]. 





Doublet 


Quartet 




14.1 MeV 42.0 MeV 


14.1 MeV 42.0 MeV 


4 


0.4662 


0.4929 


0.9794 


0.8975 




(0.4710) 


(0.4719) 


(0.9809) 


(0.8865) 


8 


0.4637 


0.4993 


0.9784 


0.9026 




(0.4670) 


(0.4985) 


(0.9794) 


(0.9050) 


12 


0.4640 


0.5014 


0.9783 


0.9030 




(0.4664) 


(0.5041) 


(0.9792) 


(0.9071) 


16 


0.4643 


0.5019 


0.9782 


0.9031 




(0.4666) 


(0.5051) 


(0.9792) 


(0.9071) 


20 


0.4644 


0.5021 


0.9782 


0.9033 




(0.4666) 


(0.5052) 


(0.9791) 


(0.9069) 


24 


0.4645 


0.5022 


0.9782 


0.9033 




(0.4666) 


(0.5055) 


(0.9790) 


(0.9071) 


28 


0.4645 


0.5022 


0.9782 


0.9033 




(0.4667) 


(0.5056) 


(0.9790) 


(0.9071) 


Ref.[8] 


0.4649 


0.5022 


0.9782 


0.9033 



A. Inelasticity and phase-shifts 

The unitarity of the 5-matrix implies that given an 
incoming channel, for instance channel 1 (n + d channel) , 
we have that l<5i«| 2 — 1; or i m other words, 

oo 

^|5 ln | 2 = l-|5 u | 2 , (24) 

n=2 

which means that an accurate calculation of the elastic 
term Sn amounts to an accurate calculation of the infi- 
nite summation of the |iSi„| 2 terms (n > 1) corresponding 
the breakup channels. 

Also, the complex value of iSn can be written in terms 
of a complex phase-shift S as: 

Sn = e m = e -^(5) e ^(S) = | <Sn | e 2iRe(5)_ (25) 

The value of |6>ii| 2 gives the probability of elastic 
neutron-deuteron scattering, and |<Sn| is what usually 
referred to as the inelasticity parameter (denoted by rj 
in 0, §)• Obviously, the closer the inelasticity to 1 the 
more elastic the reaction. In fact, for energies below the 
breakup threshold the phase-shift is real and |<Sn| = 1. 

In Table U we give the inelasticity parameter |<Sn| for 
the two laboratory neutron energies used in ||, i.e., 14.1 
MeV and 42.0 MeV. We have computed the inelasticity 
parameter for both, the doublet and the quartet states. 
We give the computed values of |5n| for different trunca- 
tions in the infinite summation (|18[) . In particular, in the 
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TABLE II: The same as Table U for Re(8). 



I - 


Doublet 


Quartet 




14.1 MeV 42.0 MeV 


14.1 MeV 42.0 MeV 


4 


105.82 


42.66 


69.04 


38.98 




(97.62) 


(28.89) 


(60.88) 


(25.29) 


8 


105.57 


41.65 


68.99 


37.95 




(99.91) 


(32.88) 


(63.24) 


(28.98) 


12 


105.53 


41.49 


68.98 


37.77 




(101.01) 


(34.80) 


(64.43) 


(31.07) 


16 


105.53 


41.46 


68.97 


37.73 




(101.68) 


(36.05) 


(65.11) 


(32.28) 


20 


105.53 


41.45 


68.96 


37.72 




(102.12) 


(36.81) 


(65.54) 


(33.03) 


24 


105.53 


41.44 


68.96 


37.71 




(102.41) 


(37.30) 


(65.84) 


(33.55) 


28 


105.53 


41.44 


68.96 


37.71 




(102.49) 


(37.90) 


(65.88) 


(33.58) 


Ref.[8] 


105.50 


41.37 


68.96 


37.71 



table we give the value K max of the asymptotic grand- 
angular quantum number associated to the last adiabatic 
potential included in the calculation (see Eq.(fl"6|)). The 
values given in the table without parenthesis have been 
obtained by using the integral relations (fll?l and (I2TJ1) . 
while the ones within parenthesis have been calculated 
from the A^ K ' and B^ K ~> matrices extracted directly from 
the asymptotic part of the three-body wave function, as 
indicated in Eq. (fT4"l) . The last row in the table gives the 
value obtained in Q. 

As seen in the table, when using the integral relations 
the agreement with the results in 0] is very good. Ac- 
tually, we obtain precisely the same result for the two 
energies in the quartet case, and a tiny difference clearly 
smaller than 0.1% in the doublet case. Furthermore, the 
pattern of convergence is rather fast, specially in the 
quartet case, for which already for if max = 8 we ob- 
tain a result that can be considered very accurate. In 
the doublet case the convergence is a bit slower, and a 
value of -Kmax of about 16 is needed. As we can see from 
the values within parenthesis, when the integral relations 
are not used, the value of |<Sn| seems to converge more 
slowly, and even if converged, the result is less accurate. 

From Eq. (l25l) we have that, together with the inelastic- 
ity parameter |iSu|, a complete specification of the matrix 
element S\± requires knowledge of the real part of the 
phase-shift Re(<5). The corresponding computed values 
are shown in Table UH where the meaning of the differ- 
ent columns is the same as in Table |U The behavior of 
Re(<5) is similar to what shown in Table fl] for |<Sn|. For 
the quartet case the same results as in Q are obtained, 
while for the doublet again a very small difference smaller 
than 0.1% is again found. The pattern of convergence is 
also similar. A K max value of 12 is already enough to get 
a quite accurate value. The main difference compared 



to the results for |<Sn| shown in Table Q] is that now the 
values of Re(<5) obtained without use of the integral re- 
lations are much less converged and much less accurate. 
The difference with the true result can reach up to 10%. 
This behavior was already observed in [l| for the phase- 
shift in a reaction below the breakup threshold. This is 
due to the fact that the hyperradius p and the Jacobi 
coordinate y entering in the asymptotic forms (I12[) and 
(|13p are equivalent only at infinity as commented at the 
end of section Hi Bl . This means that a correct extraction 
of the phase-shift from the asymptotic part of the wave 
function requires to impose the boundary condition at 
infinity, for which also infinite adiabatic terms would be 
needed. 

Therefore, for energies above the breakup threshold, 
the conclusion from Tables U and |TT] is similar to the one 
reached in [lj for 1 + 2 elastic processes, that is, the use of 
the integral relations is crucial from two different points 
of view. First, it accelerates drastically the convergence 
of the values of the 5-matrix elements, and second, they 
arc needed in order to obtain the correct result. 



IV. SOFT-CORE 4 He- 4 He POTENTIAL AND 
THE 4 He- 4 He 2 REACTION 

Once the integral relations have been proved to be 
efficient in order to describe 1 + 2 reactions above the 
breakup threshold, in this section we shall use them to 
study the atomic 4 He- 4 He2 process. In particular, we 
shall focus on three-body states with spin and parity + , 
and the possibility of using simple two-body soft-core po- 
tentials will be investigated. 

The 4 He- 4 He molecule is known to be one of the biggest 
diatomic molecules. Its binding energy has been esti- 
mated to be around 1 mK with a scattering length a 
around 190 a.u (l2l,[l3j. Different accurate investigations 
of the 4 He- 4 He interaction are available in the literature 
[T3-LT3] ■ AU of them present the common feature of a 
sharp repulsion below an interparticle distance of approx- 
imately 5 a.u.. The presence of the repulsive core is the 
source of a series of important technical difficulties. For 
instance, the wave function in the inner regions, which 
is very small due to the large potential repulsion, is de- 
cisive for the energy of the bound states or the asymp- 
totic properties of the continuum states. Therefore, the 
wave function must be calculated with high accuracy in 
this region, which typically requires a very important in- 
crease of the basis size. Furthermore, in case of using 
the adiabatic expansion method, the angular eigenvalues 
\ n (p) (see Eq.([3])) also diverge for small p, and this di- 
vergence provokes very frequent crossings between them 
that sometimes are not easy to handle 18]. 

In Ref.[9[ the possibility of using a soft-core potential 
was investigated in the context of 4 He- 4 He2 collisions be- 
low the threshold for breakup of the dimer. In particular 
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the gaussian potential suggested in [18| 
V 2b {r) = -1.227e-'- 2 / 10 - 032 



(26) 



was used (the strength of the potential is in K and the 
range in a.u.). 

The dimer properties obtained with this potential are 
given in the second column in the upper part of table Hill 
In particular, the dimer binding energy (-E2&)) the two- 
body scattering length (a), the effective range (ro), and 
the interatomic distance ((r)) are given. The two param- 
eters in the soft-core gaussian potential (|26)) were fitted 
to reproduce the scattering length and the effective range 
of the hard-core potential LM2M2 [1J|. As seen in the 
third column of Tabic [TTTI (upper part), when this is done 
the binding energy E 2 b and the interatomic distance (r) 
are also well reproduced. 

However, even if both potentials have the same two- 
body properties, when moving to the three-body states 
important differences appear. The soft-core potential 
overbinds the two bound states in 4 He3 and the atom- 
dimer scattering length is clearly smaller (second and 
fourth columns in the lower part of Table IIII[) . In fact, 
as shown in [9j , the phase-shifts obtained with these two 
potentials for the 4 He — 4 He2 elastic scattering (below the 
breakup threshold) clearly differ from each other. As an 
example, for an incident energy of 1 mK, the phase-shifts 
obtained with the gaussian and the LM2M2 potentials 
are —56 and —63 degrees, respectively. 

As also shown in Q, this anomalous behavior of the 
soft-core potential at the three-body level can be cor- 
rected by using a short-range effective three-body force, 
depending only on the hyperradius, that is added to the 
effective potential In particular we choose the simple 
three-body force 



W(p) = W e- p2/p °, 



(27) 



where the strength Wo is adjusted to reproduce the 
trimer ground-state binding energy obtained with the 
LM2M2 potential. When this is done the results are 
fairly independent of the range parameter po, at least 
within a reasonable value from po = 4 a.u. to po = 10 
a.u.. The results given in the third column (lower part) 
of Table IIIII have been obtained after inclusion of the 
three-body force. The precise values of Wq and po used 
are given in the first row of the lower part of the table. 
As we can see, once the binding energy of the ground 
state has been corrected, the binding energy of the ex- 
cited state and the atom-dimer scattering length auto- 
matically agree with the corresponding values obtained 
with the LM2M2 potential. Furthermore, as seen in Fig. 3 
of @ , the low energy phase-shifts are also corrected when 
the effective three-body force is included. 

The importance of the inclusion of the three-body force 
is also seen when investigating the 4 He — 4 He2 reaction 
for incident energies (Ei = E — E 2 b) above the threshold 
for breakup of the dimer, i.e., Ei > \E 2 b\ (or E > 0). In 
Table HVl we give the inelasticity (|«Sn|) and the real part 



of the phase-shift (Re(5)) for energies E = 5 and E = 25 
mK. As in Tables Q] and HH K maK is the grand-angular 
quantum number associated to the last adiabatic term 
included in the expansion (JXJ) . The boldface numbers 
have been obtained when the three-body force has been 
included in the calculation (Wo = 18.314 K and po = 6 
a.u.), and the results within parenthesis have been ob- 
tained without the three-body force. 

As we can see, the pattern of convergence is similar 
to the one observed in Tables U and |H] for the neutron- 
deuteron reaction. A K m&K value of around 12 is enough 
to get a rather well converged inelasticity, while Re(S) 
requires a few more adiabatic terms in order to reach 
convergence. Similarly to what found in Q for ener- 
gies below the breakup threshold, the inclusion of the 
three-body force gives rise to relevant changes in the com- 
puted values. These changes are particularly noticeable 
for Re((5), which for the two energies under consideration 
increases up to 6 degrees when the three-body force is 
not included in the calculation. Also, we observe that for 
E = 5 mK the breakup probability (1 — |tSn| 2 ) is still 
rather small, clearly smaller than 1% (with three-body 
force), while for E = 25 mK this probability rises up to 
17% (12% without three-body force). 

With the iS-matrix in hand, we can now compute the 
dissociation rate for the 4 He — 4 He2 collision. The an- 
alytic form of this rate is given in Ref.|20], and for + 
states it becomes: 



where 



Mi, 23 



Tin 



Ml, 23^1, 23 



(l-|5n| 2 ), 



2M He 



k 2 

K l,23 



2^1,23-^1 

h 2 



(28) 



(29) 



Mn e is the mass of the 4 He atom, and Ei — E — E 2 b is 
the incident energy in the center of mass frame (E 2 b is 
the binding energy of the 4 He2 dimer). 

In FigJ3] we show the dissociation rate as a function 
of the three-body energy E. The thin-dashed and thin- 
solid lines are the results obtained with the gaussian soft- 
core potential (|26p without and with the additional three- 
body force (|27[) . respectively. The low energy behavior 
of the dissociation rate follows the E Km+2 rule derived 
in [2l| , where K m is the smallest grand-angular quantum 
number associated to the continuum adiabatic channels 
(K m = in our case of three indistinguishable bosons 
coupled to J n = 0+). 

From the figure we can see that the effect of the three- 
body force is quite important. In fact, for small energies 
the three-body force reduces the rate by a factor of 2. At 
higher energies, beyond 10 mK, the effect is the opposite, 
and the three-body force increases the rate by a factor 
close to 1.5. 

In Ref.[2(| the same dissociation rate has been com- 
puted. The corresponding curve is shown in Fig|3] by 
the dotted curve. As we can see, there is an important 
difference compared to our calculation. Except at high 
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TABLE III: Helium dimer and helium trimer properties obtained with the soft-core potential in Bq, (126|l . the hard-core potential 
LM2M2 in Ref.[l4j|. the modified soft-core potential (soft-core 2) in Eg. (|30|l . and the hard-core SAPT potential in Ref.[T3], 
respectively. The two-body properties are given in the upper part of the table, where E^b, a, ro, and (r) are the dimer binding 
energy, the scattering length, the effective range, and the average interatomic distance, respectively. The three-body properties 
are given in the lower part, where E^ b ' s '\ E^ xc '\ and ao are the binding energy of the trimer ground state, the binding energy 
of the trimer excited state, and the atom-dimer scattering length, respectively. For the soft-core and soft-core-2 potentials the 
column corresponding to the trimer properties has been split in two. The results given in the left part have been obtained 
without inclusion of any effective three-body force. The ones given in the right part are obtained after inclusion of a three-body 
force fitted to match the binding energy of the trimer ground state obtained with the LM2M2 [T^ and SAPT [l7],[23| potentials, 
respectively. The parameters used for the three-body potential are denoted by Wo and po (see Eq.(|27[l). All the distances are 
given in a.u. and the energies in mK, except Wo which is given in K. 







Soft-core 


LM2M2 [19] 




Soft-core 2 


SAPT [17, 201 


E 2b (mK) 




-1.296 


-1.302 




-1.554 


-1.564 


a (a.u.) 




189.95 


189.05 




174.09 


173.50 


ro (a.u.) 




13.85 


13.84 




13.80 


13.79 


(r) (a-u.) 




98.4 


98.2 




90.5 


90.3 


(Wo, po) (K, a.u.) 


(o,-) 


(18.314, 6) 




(o,-) 


(17.760, 6) 




4r° (mK) 


-150.0 


-126.4 


-126.4 


-154.9 


-130.9 


-130.9 


sir ( mK ) 


-2.467 


-2.287 


-2.265 


-2.805 


-2.612 


-2.588 


a (a.u.) 


165.9 


210.6 


224.3 


181.7 


226.0 


226.8 



TABLE IV: Inelasticity (|iSn|) and real part of the phase- 
shift (Re(<5)) for the 4 He- 4 He2 collision at three-body energies 
(above threshold) E — 5 mK and E = 25 mK. The value of 
K m &x is the AT- value associated to the last adiabatic potential 
included in the calculation. The numbers within parenthesis 
have been obtained without inclusion of the three-body force. 



-^max 


|Sn| 


Re(<5) 




E=5mKfi 


= 25 mK 


E = 5mK E 


= 25 mK 


4 


0.9988 


0.9351 


69.30 


35.09 




(0.9946) 


(0.9645) 


(75.52) 


(40.59) 


8 


0.9988 


0.9111 


69.23 


34.80 




(0.9946) 


(0.9403) 


(75.45) 


(40.28) 


12 


0.9989 


0.9104 


69.20 


34.63 




(0.9947) 


(0.9394) 


(75.41) 


(40.09) 


16 


0.9989 


0.9110 


69.17 


34.60 




(0.9947) 


(0.9402) 


(75.39) 


(40.05) 


20 


0.9989 


0.9110 


69.16 


34.59 




(0.9947) 


(0.9402) 


(75.38) 


(40.04) 


24 


0.9989 


0.9109 


69.15 


34.58 




(0.9947) 


(0.9402) 


(75.38) 


(40.04) 


28 


0.9989 


0.9109 


69.15 


34.58 




(0.9947) 


(0.9402) 


(75.37) 


(40.04) 


40 


0.9989 


0.9109 


69.15 


34.58 




(0.9947) 


(0.9402) 


(75.37) 


(40.04) 



energies (E > 10 mK), where our result (with three-body 
force) and the one in [20J basically coincide, for small en- 
ergies our rate is about a factor of 3 bigger. 

However, we have to note that the hard-core two-body 
potential used in [2(| gives rise to somewhat different 
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FIG. 3: (Color online) Dissociation rate ( - Eq. (l28l) ') for the 
4 He + 4 He 2 -> 4 He + 4 He + 4 He reaction as a function of 
the three-body energy E. The thin and thick curves have 
been obtained with the soft-core potential (|26|l and (|30|l . re- 
spectively (see Table iHTjl . The dashed and solid curves are the 
results without and with inclusion of an effective three-body 
force. The dotted curve is the result given in (2pj ] . 



dimer properties compared to the soft-core potential (T2rj|) 
and the LM2M2 potential Q. In [Io[ they have used the 
potential based on the symmetry- adapted perturbation 
theory (SAPT) derived in 



17| . The two-body properties 



obtained with this potential are given in the last column 
(upper part) of Table IIIII The two-body dimer is 20% 
more bound than with the LM2M2 potential, and there- 
fore the scattering length and the interatomic distance 
are smaller than in the LM2M2 case. 

To investigate the sensitivity of the dissociation rate 
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FIG. 4: (Color online) The same as Fig|3] for the recombi- 
nation rate (Eq. (|31[) ') for the process 4 He + 4 He + 4 He — > 
4 He 2 + 4 He. 



to the details of the two-body interaction we have then 
constructed a second gaussian soft-core potential repro- 
ducing the two-body properties of the SAPT potential. 
This potential takes the form: 



V 2b (r) = -1.234e- r2 / 10 - 032 , 



(30) 



where the strength is given in K and the range in a.u.. 
The corresponding two-body properties are given in the 
fourth column (upper part) of Table Mil under the label 
"soft-core 2" . Again, when moving to the three-body sys- 
tem, this new soft-core potential presents the same defi- 
ciencies as the previous one, i.e., overbinding of the three- 
body bound states and a too small atom-dimer scattering 
length. As before, this problem is solved by inclusion of 
the three-body force, whose parameters are again fitted 
to reproduce the ground state binding energy of the he- 
lium trimer provided by the SAPT potential. The three- 
body properties with and without three-body force, as 
well as the parameters used for the three-body force are 
given in the last three columns (lower part) of Table Hill 

Making use of this new soft-core potential and the cor- 
responding three-body force, we can compute again the 
dissociation rate ([28]) . The results are shown in Figj3] 
by the thick-dashed curve (no three-body force) and the 
thick-solid curve (with three-body force). The new rates 
are now smaller than the ones obtained with the previous 
soft-core potential. Furthermore, when the three-body 
force is included, and therefore not only the two-body 
properties but also the three-body ones agree with the 
ones obtained with the SAPT potential, our dissociation 
rate and the one given in [2(| agree well for the whole 
range of energies. Only a small difference can be seen at 
about 1 mK. We can therefore see that the dissociation 
rate is very sensitive to the details of the two-body inter- 
action. An additional 20% binding of the dimer, and the 
corresponding decrease in the scattering length, leads to 
a factor of 3 decrease in the rate. 

The same conclusions are reached when investigating 



the recombination rate corresponding to the inverse pro- 
cess 4 He 4- 4 He + 4 He — s> 4 He2 + 4 He. This rate is given by 
the same iS-matrix elements as in the dissociation case, 
and its analytical form for three identical bosons with 
angular momentum and parity + is given by [20|: 



= 3! 



32/lTT 2 



(1-15 



ii 



where 



ML 



k 2 = 



(31) 



(32) 



Mue is the mass of the 4 He atom, and E is the total 
three-body energy. 

The computed recombination rate is shown in Figj4] 
The meaning of the curves is as in Fig(3J In Ref.[2l| the 
recombination rate K% was proved to behave at low ener- 
gies as E Km {K m is the smallest grand-angular quantum 
number associated to the continuum adiabatic channels), 
which means that in our case, where K m = 0, K3 should 
be constant at very low energies, as observed in the fig- 
ure. Again, inclusion of the three-body force reduces the 
rate by a factor of 2 to 3 at small energies (compare the 
dashed curves and the corresponding solid curves in the 
figure). Also, when the soft-core potential (|3"0)l together 
with the three-body force is used, the computed recombi- 
nation rate agrees well with the rate given in [2(| (dotted 
curve in the figure). This fact reveals again the impor- 
tance of the fine details of the two-body interaction. 

Finally, it is important to emphasize that in this work 
the iS-matrix has been obtained by use of the integral re- 
lations dl~9|) and ([20]) . As explained, they have the great 
advantage of needing only the internal part of the wave 
functions. In fact, the integrals involved in the calcula- 
tions shown in this section can be safely performed with a 
maximum value for the hyperradius p of 4000 a.u.. Even 
less could be enough. However, in [20(, where the adia- 
batic expansion was also used, the radial wave functions 
in Eq.((T|) had to be expanded up to 5 x 10 5 a.u., which is 
definitely very delicate from the numerical point of view. 



V. SUMMARY AND CONCLUSIONS 

In this paper we have extended the use of the integral 
relations derived in 0,0 to describe 1+2 reactions above 
the threshold for breakup of the dimer. As in 0, 0, the 
integral relations are used in combination with the adi- 
abatic expansion method in order to construct the trial 
wave function. The main consequence of moving up to 
energies above the breakup threshold is that the 5-matrix 
describing the process has now infinite dimension. This 
was not the case for energies below the threshold, where 
the elastic, inelastic, or transfer channels were described 
by a finite (and small) number of adiabatic terms. 

The applicability of the method has been tested with 
the neutron-deuteron reaction, for which benchmark cal- 
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culations are available. The agreement with these calcu- 
lations is good, and the pattern of convergence is sim- 
ilar to the one found in [3, 0] for energies below the 
breakup threshold. The integral relations accelerate the 
convergence significantly, and typically about 10 adia- 
batic terms are enough to get a rather well converged 
result. 

The method has then been used to investigate reac- 
tions involving three helium atoms. In particular we have 
focused on the possibility of using soft-core atom-atom 
potentials in order to describe the full process. These 
potentials permit to avoid all the technical problems aris- 
ing from the use of more sophisticated potentials where 
a hard-core repulsion is always present. However, as al- 
ready shown in |9[ , even if the soft-core potentials repro- 
duce properly the dimer properties, an effective three- 
body force is needed in order to reproduce as well the 
properties of the trimcr bound states and the atom-dimer 
scattering length obtained with the hard-core potentials. 

We have found that the three-body force also modifies 
significantly the inelasticity and the phase-shift for en- 
ergies above the breakup threshold. In fact, when com- 
puting the reaction rates for dissociation of the dimer 
( 4 He + 4 He2 —5- 4 Hc + 4 Hc + 4 He) and for the recom- 
bination process ( 4 He + 4 He + 4 He — > 4 He2 + 4 He) the 
three-body force decreases the rates even by a factor of 
3 at low energies. 

These two rates are also very sensitive to the details of 
the two-body interaction. We have found that a two- 
body potential, like the SAPT potential, providing a 



dimer state about 20% more bound than the one ob- 
tained with the LM2M2 potential, also reduces the rates 
by a factor of around 3 at small energies. The soft- 
core potential reproducing the two-body properties of the 
SAPT potential, together with the corresponding three- 
body forced designed to reproduce the trimer properties 
as well, gives then rise to reaction rates in very good 
agreement with the ones of the SAPT potential. 

In summary, the integral relations are also useful in 
order to describe reactions above the breakup thresh- 
old. They accelerate significantly the convergence of the 
5-matrix terms. For reactions involving hard-core two- 
body potentials, the use of soft-core potentials with the 
same two-body properties are a very good alternative, 
provided that they are used together with an effective 
three-body force designed to fit as well the bound state 
three-body energies. When this is done, the reaction 
rates obtained with the hard-core and the soft-core po- 
tentials agree pretty well. These rates are very sensitive 
to the details of the two-body interaction. Small varia- 
tions of the dimer properties can produce sizable changes 
in the reaction rates at low energies. 



Acknowledgments 

This work was partly supported by funds provided by 
DGI of MINECO (Spain) under contract No. FIS2011- 
23565. 



[1] P. Barletta, C. Romero- Redondo, A. Kievsky, M. Viviani, 

E. Garrido, Phys. Rev. Lett 103, 090402 (2009). 
[2] C. Romero- Redondo, E. Garrido, P. Barletta, A. Kievsky, 

M. Viviani, Phys. Rev. A 83, 022705 (2011). 
[3] F.E. Harris, Phys. Rev. Lett. 19, 173 (1967). 
[4] A.R. Holt and B. Santoso, J. Phys. B 5, 497 (1972). 
[5] A. Kievsky, M. Viviani, P. Barletta, C. Romero- Redondo, 

E. Garrido, Phys. Rev. C 81, 034002 (2010). 
[6] P. Barletta and A. Kievsky, Few-Body Syst. 45, 25 

(2009). 

[7] J.L. Friar el al, Phys. Rev. C 42, 1838 (1990). 

[8] J.L. Friar, GL. Payne, W. Glockle, D. Hiiber, H. Witala, 

Phys. Rev. C 51, 2356 (1995). 
[9] A. Kievsky, E. Garrido, C. Romero-Redondo, P. Barletta, 
Few-body Syst. 51, 259 (2011). 
[10] F.M. Pen'kov, J. Exp. Theor. Phys. 97, 485 (2003). 
[11] E. Nielsen, D.V. Fedorov, A.S. Jensen, and E. Garrido, 

Phys. Rep. 347, 373 (2001). 
[12] F. Luo, C.F. Giese, and W.R. Gentry, J. Chem. Phys. 
104, 1151 (1996). 



[13] R.E. Grisenti, W. Schollkopf, J.P Toennies, GC. 
Hegerfeldt, T. Kohler, and M. Stoll Phys. Rev. Lett. 85, 
2284 (2000). 

[14] Aziz R.A. and Slaman M.J., J. Chem. Phys. 94, 8047 
(1991). 

[15] K.T. Tang, J.P. Toennies, and C.L. Yiu, Phys. Rev. Lett. 

74, 1546 (1995). 
[16] T. Korona, H.L. Williams, R. Bukowski, B. Jeziorski, and 

K. Szalewicz, J. Chem. Phys. 106, 5109 (1997). 
[17] M. Jeziorska, w. Cencek, K. Patkowsky, B. Jeziorski, K. 

Szalewicz, J. Chem. Phys. 127, 124303 (2007). 
[18] E. Nielsen, D. V. Fedorov and A.S. Jensen, J. Phys. B 

31, 4085 (1998). 
[19] E.A. Kolganova, A.K. Motovilov, and W. Sandhas, Phys. 

Part, and Nucl. 40, 206 (2009). 
[20] H. Suno and B.D. Esry, Phys. Rev. A 78, 062701 (2008). 
[21] B.D. Esry, C.H. Greene, and H. Suno, Phys. Rev. A 65, 

010705 (2001). 



